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^^ Abstract 

Pm ' The interplay between space and evolution is an important issue in population 

r^^ \ dynamics, that is in particular crucial in the emergence of polymorphism and spa- 

f~| ' tial patterns. Recently, biological studies suggest that invasion and evolution are 

^ I closely related. Here we model the interplay between space and evolution starting 

with an individual-based approach and show the important role of parameter scalings 
on clustering and invasion. We consider a stochastic discrete model with birth, death, 
competition, mutation and spatial diffusion, where all the parameters may depend 
both on the position and on the trait of individuals. The spatial motion is driven 
by a reflected diffusion in a bounded domain. The interaction is modelled as a trait 
^b I competition between individuals within a given spatial interaction range. First, we 

(^ • give an algorithmic construction of the process. Next, we obtain large population 

^^ I approximations, as weak solutions of nonlinear reaction-diffusion equations with Neu- 

^^ • mann's boundary conditions. As the spatial interaction range is fixed, the nonlinearity 

is nonlocal. Then, we make the interaction range decrease to zero and prove the con- 
vergence to spatially localized nonlinear reaction-diffusion equations, with Neumann's 
boundary conditions. Finally, simulations based on the microscopic individual-based 
C^ ■ model are given, illustrating the strong effects of the spatial interaction range on the 

emergence of spatial and phenotypic diversity (clustering and polymorphism) and on 
the interplay between invasion and evolution. The simulations focus on the qualitative 
differences between local and nonlocal interactions. 
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1 Introduction 

The spatial aspect of population dynamics is an important ecological issue that has 
been extensively studied (Murray ,31], Durrett and Levin 12:, Tilman and Kareiva |44j . 
McGlade j26j . Dieckmann et al. [9J- It is in particular crucial in environmental prob- 
lems, such as spatial invasions and epidemics (Mollison ^], Murray [HJ, Rand et al. 
j36j . Tilman and Kareiva 44,, Lewis and Pacala US]), and clustering or agglomeration of 
the population, i.e. the organization as isolated patches (Hassel and May ^\, Hassel and 
Pacala [10], Niwa [32], Flierl et al. ^5, Young et al. 05]). The combination of space and 
phenotype is also known for a long time to have important effects on population dynamics 
(Mayr [27], Endler ^S])- In particular, it can strongly favor the coexistence of several types 
of individuals and the emergence and stability of polymorphism (Durrett and Levin |llj . 
Dieckmann and Doebeli [7|). More recently, several biological studies (Thomas et al. |43| . 
Phillips et al. i33j) observed that classical models could underestimate the invasion speed 
and suggested that evolution and invasion are closely related. Namely, the evolution of 
morphology can have strong impact on the expansion of invading species, such as insect 
species ([IS]) or cane toads (|S3])- In this context, the study of space-related traits, such 
as dispersal speed (Prvost [SI], Desvillettes et al. 0), or sensibility to heterogeneously 
distributed resources (Bolker and Pacala [S], Grant and Grant 18 ), is fundamental. 

In this paper, we propose and construct stochastic and deterministic population models 
describing the interplay between evolution and spatial structure. We show how helpful 
these models can reveal to understand and predict several specific behaviors concerning 
clustering and invasion. 

We study the dynamics of a spatially structured asexual population with adaptive 
evolution, in which individuals can move, reproduce with possible phenotypic mutation, 
or die of natural death or from the competition between individuals. The individuals are 
characterized both by their position and by one or several phenotypical adaptive traits 
(such as body size, rate of food intake, age at maturity or dispersal speed). The interaction 
is modelled as a trait competition between individuals in some spatial range. Our approach 
is based on a stochastic microscopic description of these individuals' mechanisms, involving 
both space and traits. This approach has already been developed in simpler ecological 
contexts. For populations undergoing dispersal, Bolker and Pacala [U [S] and Dieckmann 
and Law [H], offered the first microscopic heuristics and simulations. Their individual- 
based model has been rigourously constructed in Fournier and Meleard ^H] • If one thinks of 
the dispersion in the physical space as a mutation in a trait space, this model translates into 
an evolutionary model. The generalization to adaptive population with general mutation 



and competition phenomena is achieved by Champagnat, Ferriere and Meleard OHj. In 
these papers, different large population deterministic or stochastic approximations have 
been obtained, depending on several scalings on the microscopic parameters. 

The basic mechanisms of the population dynamics we consider combine spatial motion 
and evolutionary dynamics (Sectionj^l). The birth, mutation and death parameters of each 
individual depend on its position and trait. An offspring, appearing at the same position 
as its progenitor, usually inherits the trait value of the latter, except when a mutation 
causes the offspring to take an instantaneous mutation step at birth to new trait values. 
As soon as it is alive, an individual moves in the spatial domain according to a reflected 
diffusion process. Moreover, each individual dies because of natural death or is eliminated 
in the competition (selecting the fittest traits) between individuals living in a given spatial 
range 6 > 0. 

Section 121 starts with the algorithmic construction of a stochastic Markov point process 
whose generator captures the individual migration and ecological mechanisms in the pop- 
ulation. Then the existence of this measure-valued process and its martingale properties 
are proved under some moment condition on the initial data. 

Next (Section 0J), we study approximations of this model based on large-population 
limits. We consider a large number N of individuals at initial time and assume that a 
fixed amount of available resources has to be partitioned between individuals. When N 
tends to infinity, the conveniently normalized point process converges to a deterministic 
finite measure, solution of a nonlinear nonlocal integro-differential equation with Neu- 
mann's boundary conditions, parameterized by the spatial range. The proof is based on 
the martingale properties of the process and on limit theorems for measure-valued jump 
processes. We moreover prove that for sufficiently smooth and non-degenerate diffusion 
coefficients, assuming that the initial condition has a density, the limiting measure has at 
each time a density with respect to the Lebesgue measure. That is due to the regularizing 
effect of the refiected diffusion process. The proof mainly uses analytic tools, and is based 
on the mild formulation of the limiting nonlinear equation. 

In Section |SJ we study the behavior of this density function as the interaction range 
tends to 0. We show its convergence to the solution of a spatially local nonlinear integro- 
differential equation with Neumann's boundary conditions. This equation has been in- 
troduced and studied in Prevost [35 in an analytic point of view (see also Desvillettes et 
al. [HI). In this spatially local case, numerical simulations by finite element methods are 
given and show the infiuence of diffusion and mutation parameters on the invasion of the 
domain by the population. 



In Section IHJ we give simulations of the microscopic process illustrating the time- 
dependent interplay between space and adaptation. We address the effect of the popula- 
tion size, and the crucial role of the interaction range with respect to spatial organization 
(clustering) and polymorphism. We focus on the qualitative differences between nonlocal 
and local interactions. In a first example, we show that, when migrations and mutations 
are not too strong, a large interaction range induces a spatial organization of the popu- 
lation as a finite set of isolated clusters, as assumed in classical metapopulation models 
(|llj). Such a spatial organization is related to the ecological notion of "niches" (different 
types of individuals settle different regions of space, Roughgarden jSH])- Conversely, for 
sufficiently small interaction range, the clustering phenomenon is no more observed. Next, 
we propose another example where a similar phase transition occurs for spatial cluster- 
ing and in which the critical interaction range can be identified. In our last example, 
we investigate a model describing the invasion of a species with evolving dispersal speed 
(as in i6j). The diffusion coefficient and the trait are assumed to be proportional and a 
triangular invasion pattern is observed, indicating that the invasion front is composed of 
faster individuals ([SS])- 

Notation 

The individuals live in the closure of a bounded domain X of W^ of class C^ and their 
trait values belong to a compact set U of M'^. 

- For X E dX, we denote by n{x) the outward normal to the boundary dX at point x. 

- For a sufficiently smooth function / and (x,u) G dX x U, we denote by dnf{x,u) the 
scalar product Vxf{x,u) ■ n{x). 

- We denote by Cq' the space of measurable functions /(x, u) of class C^ in x and bounded 
in u satisfying dnf{x, u) = for all (x, u) £ dX x U and by Cq' the subspace of functions 
f{x,u) which are moreover continuous in u. 

- For each p > 1, the L^-norm on X xU is denoted by || • ||p. 

- We denote by Mp{X x U) the set of finite measures on X x U, endowed by the weak 
topology, and by M the subset of Mf{X x U) composed of all finite point measures, that 
is 

M = I ^(^(xi,n*)i "- G N, x^,...,x^ € X , u\...,u" ^u\ 

where 6tx,u) denotes the Dirac measure at {x,u). (If n = 0, one obtains by extension the 
null measure). For any z^ G Mf{X x U) and for any measurable function f on X x U, we 
write indifferently {u, f) or /^^^ fdu. li u = X]"^^ (^(^»,n»), then {u, f) = X^ILi f{x\u') . 

- We denote by B([0, oo),Mi7'(^ x U)) the Skorohod space of left limited and right con- 



tinuous functions from M+ to Mp[X x U), endowed with the Skorohod topology. 
- The constant C will be a constant which can change from line to line. 

Remark 1.1 Let us remark that the space of C'^{X) -functions with a vanishing normal 
derivative is dense, for the uniform norm, in C{X). Indeed, let us consider the Cauchy 
problem for the parabolic differential equation ^{t,x) = Au{t,x) ; t > ; x G A" with the 
boundary condition ^(t, x) = ; t > ; x G dX . Since X is of class C^ , we may apply 
Sato-Ueno HW Theorem 2.1. There exists a smooth fundamental solution q{t,x,y) to this 
problem and each f € C{X) is the uniform limit of the sequence J-j^q{t,x,y)f(y)dy of 
C'^ (X) -functions with vanishing normal derivative, as t tends to 0. 

We easily extend this result and show that the space Cg' is dense in the space of 
continuous functions on X x U. 

2 The model 

Let us now describe the evolutionary process we are interested in. The population will be 
described at any time by a finite point measure z^ S A^. Each individual, characterized by 
its position and trait {x,u), may move, give birth or die, as described below. 

1. The migration is described as a diffusion process normally reflected at the boundary 
of the domain X. Biologists usually assume that the random behavior is isotropic, 
so the diffusion matrix is chosen with the form m{x, u)ld (Id is the identity matrix 
on M ) and the nonnegative coefficient m{x,u) (depending on the position x and 
the trait value u), is the diffusion coefficient. We moreover model the environment 
heterogeneity (resources, topography, external effects,. . . ) by a drift term driven by 
a M -vector b{x,u). 

2. Births and mutations. We consider a population with asexual reproduction. An 
individual with position x and trait u can give birth either to a clonal child at rate 
A(x,u), or to a mutant with trait v at rate M(x,u,v), both at position x. 

3. The death rate /i of an individual depends on its position x and trait u and on the 
spatial and phenotypic competition with the individuals located around x. Let us 
call (5 > the range of this spatial interaction. 



For a population v = Yl^=i \x\u^) ^ -^i the death rate is given by 

^{x, u,I W -k v{x, u)) = fiix,u, I [x — y)W{u — v)u{dy, dv) 

\ JxxU 



^ (x, u, ^iHx- x')W{u -u')j. 



The function fi{x, u, r) is assumed to be measurable on A' x ^ x M. 

This interaction assumes that spatial and phenotypic interactions are independent, 
which is realistic in many biological situations. One could of course consider a more com- 
plicated interaction. Since our ultimate goal is to make the spatial interaction range go 
to zero, we have chosen this particular form. 

Hypotheses (H): 

1) The coefficients m{x, u) and b(x, u) depend Lipschitz continuously on the position and 
measurably on the trait, and there exist constants m* > and 6* > such that for all 
{x,u) £ X xU 

< m(x,u) < m* 

(2.1) 
\bix,u)\ <b*. 

2) It is natural from a biological point of view to assume that all birth rates are bounded. 
There exists A* such that 

< A(x, u) < A*, V(x, u) £X xU. (2.2) 

The kernel M is nonnegative and symmetric in (u, v) for each x £ X and 

sup M{x,u,v) = M*{v) G L\U). (2.3) 

3) There exists a positive constant /i* such that 

V(x, u,r) eX xU xR, < fi{x, u, r) < /i*(l + \r\) (2.4) 

4) For each 5 > 0, the spatial kernel I is nonnegative and bounded and for each x €z X , 

[ I^{x - y)dy = 1. 
Jx 

5) The competition kernel W is nonnegative and bounded on W . 

Let us remark that if / is proportional to l{|x|<5}) then (H-4) means that the interac- 
tion is proportional to the surface in X around x. This is a natural biological assumption. 
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especially if x lies on the boundary of X. We will later assume that the measure I {y)dy 
weakly converges to the Dirac measure 60 as S tends to 0. 

Hypotheses (H) will be assumed in all the sequel. They imply in particular that for 
each u G Mf{X x U) and each (x, u) G X xU, 

fi{x, u, I^W • v{x, u)) < /u*(l + ||/VF||oo(j^, 1)) (2.5) 

which yields 

lJi{x,u,I W ■ki'{x^u)) + \{x,u) + I M{x,u,v)dv 

Ju 

< fi*{l + \\I^W\\oo{iy, 1)) + A* + ||M*||i < Cs{{u, 1) + 1) (2.6) 

and the total jump rate for a population i^ is bounded by 

CsW,l)i{i^,l) + l). (2.7) 

We are interested in the evolution of the stochastic point process (t'i), taking its values 
in Ai and describing the evolution of the population at each time t. We define 

Nt 

Nt ^ N standing for the number of living individuals at time t, Xl,...,X^ ' describing 
their positions (in X) and U}, ...,11^ * their trait values (in U). 

The dynamics of the population can be roughly summarized as follows. The initial 
population is characterized by a measure vq ^ ^A at time t = 0, and any individual located 
at X G ^ with trait u at time t has four independent exponential clocks: a "clonal re- 
production" clock with parameter A(x,u), a "mutant reproduction" clock with parameter 
M{x,u,v), and a "mortality" clock with parameter /i(a;,u, ^^^ /(x — X/)M^(u— t^^')). If 
the "mortality" clock of an individual rings, then this individual disappears; if the "clonal 
reproduction" clock of an individual rings, then it produces at the same location an indi- 
vidual with the same trait as itself; if the "mutant reproduction" clock of an individual 
rings, then it produces at the same location an individual with characteristics {x^v). 

The living individuals evolve in the domain, according to diffusion processes with 
diffusion coefficient m{x,u) and drift b{x,u), normally reflected at the boundary of X. 

The measure- valued process (t't)t>o is a Markov process whose infinitesimal generator 
L captures this dynamics. This generator is the sum of a jump part Li corresponding 
to the phenotypic evolution and of a diffusion part L2. The generator Li is defined for 



bounded and measurable functions (j) from M. into M and for u = Yl^=i ^(x\u^) by 

i=l 

+ / ^{(pii^ + S^^^^^)) - (t>{iy))M{x\u\v)dv 

+ ^((/>(i/ - <5(,.,„.)) - ct>{u))ix{x\u\ I'W * v{x\ u')). (2.8) 

A standard class of cylindrical functions generating the set of bounded and measurable 
functions from Mf{X x lA) into R is the class of functions 

Ff{u)=F{{vJ)), (2.9) 

for bounded and measurable functions F and /. 

For such functions Fj, with F G C'^(R) and f £ Cq' , the diffusive part L2 of the 
generator can easily be deduced from Ito's formula. Its form is similar to the one obtained 
in the whole space for branching diffusing processes (cf. Roelly-Rouault [HSI) and is given 

by 

L2Ff{u) = {u, mA,./ + b.V.J)F'{{u, /)) + {u, m\Vjf)F"{{u, /)). (2.10) 

Hence, 

LFf{i^) = LiFf[u) + L2Ff(u) 

\x{x,u){F{{,.,f) + f{x,u))-F{{u,f))) 

+ / {F{{vJ) + f{x,v)) - F{{v, f)))M{x,u,v)dv 
Ju 

+ fx{x,u,I^W *u{x,u)){F{{uJ) - f{x,u)) - F{{uJ))) 
+ (m(x, u)A^f{x, u) + b{x, u).V^f{x, u))F'{{ij, /))) 

+ mix,u)\V,f{x,u)\^F"{{u,f))] v{dx,du) (2.11) 



3 Construction of the particle system and martingale prop- 
erties 

In this section, we construct a Markov process on the path space Ii{[Q, 00), Mp{X x li)) 
with infinitesimal generator L. Then we prove some martingale properties satisfied by this 



process, which are the key point to obtain large population approximations. 

Let us firstly present an iterative construction of the process, which gives an effective 
simulation algorithm, if combined with a diffusion simulation step such as an Euler scheme 
for reflected diffusions (see Lepingle j21], Gobet [ITj and Sectional). 

The initial number of individuals is equal to some natural integer N gN* and the vector 
of random variables (Xq, Uq) = {Xq, C/o)i<i<Af G (^ x U)^ denotes the position and trait 
values of these individuals. More generally, we denote by Nt the number of individuals 
at time t and by (Xt, Ut) the vector of their positions and traits. Let us introduce the 
following sequences of independent random variables, independent of {Xq, Uq). 

- {B^' )kjeN* are d-dimensional Brownian motions, 

- {9k)k are uniform random variables on [0, 1], 

- (Vfc)fc take values in U with law TrjjW-dv, 

- {Tk)k are exponential random variables with law Cse~^^lt>o. (The constant Cs is 
defined in (HSJ). 

The system is obtained inductively for fc > 1 as described below. We set Tq = and 
A^'o = A^. Assume that {Tk-i.,Nk-i-,XTj,_^^-,UTf._^) are given. If A^fc_i = 0, then i^t = for 
all t > Tk-i- If not, let 

• Tk = Tk-i + j^ — 7^^ — xjy. Notice that -^ — tj^ — -j-pr represents the time between 
possible jumps for Nk-i individuals and that CsNk-i{Nk-i+l) gives an upper-bound 
on the total jump rate for a population with Nk~i individuals, as seen in 1)2. 7(1 . 

• On the time-interval [Tk-i,Tk), the number of particles remains equal to Nk~i, 
their trait values are equal to U^ , 1 < i < Xk-i and their positions (X^, 1 < j < 
Nk-i) evolve according to the following stochastic differential equation with normal 
reflection : Vt G [Tk-i,Tk], 



X} e X, 

ft ft 


-kl. 




XI = Xi,^^^ + / ^2m{XlU^^JdB^/^ + / h{Xl Ul^jds - 
W\t = f l^^^^^g^^dWU ;H= f n{Xi)dW\s. 


) 

(3.1) 



At time Tk, one chooses at random an individual Ik = i uniformly among the 
Nk-i individuals living during the time-interval [Tk^i,Tk)- Its position and trait 
aieiXJr^^U^^J. 



-It [) < 0k < Q(JVfc„i+l) - ^l(^Tfc,t/Tfe_J, 

then the individual i dies and A'^fc = Nk-i — 1- 

- If 01(Xt„C/t,_J < 0. < e\{XT,.UT,_,) + c'CnZXi) = ^2(^T„^T,_J, then 
the individual i gives birth to an offspring with characteristics {X^ , U^ ) 

and Nk = Nk-i + 1. 

M(X^ ,U^ _ ,Vfc)||M*|li 

- If6'^(XTfc,C/Tfc_J <ek< 0|(Xr^,C/Tfc_i)H M*(Vfc)4(JVfc-i+l) " ^3(^Tfe,^Tfc_i, Vfe) 

then the individual i gives birth to a mutant offspring with trait Vk at the po- 
sition XL and N^ = Nk-i + 1. 

- If 6lfc > 6'|(XTfe, ?7Tfe_i, Vfe), nothing happens and Nk = iV^-i- 



The total number A''^ of individuals at time t is equal to Nt = J2k>o ^{Tk<t<Tk+i}^k, 



and ut = Y.k>Q hn<t<n+^} Er=\ ^(Xi,w^^) = E£i ^. 



This stochastic individual-based process v can be rigorously expressed as solution of 
a stochastic differential equation driven by d-dimensional Brownian motions (i?*)igN* and 
the M-i- X N X [0, 1] x ^-valued multivariate point process 

Q{dt,di,d9,dv) = ^5(Tkj^^e^y^){dt,di,de ,dv) 

k>l 

associated with the birth, mutation and death of individuals. We will prove its existence 
on M+, deduced from moment properties, and develop some martingale properties that we 
will use below. 

Let us consider fo G M.. For each Cq' -function /, we define the process {vt, f) as 
solution of the stochastic differential equation 

(t't,/) = (z^o,/) + / {i^r,m{x,u)Axf + b{x,u).Vxf)dr 
Jo 

„tK-,i> . 

+ / V J2m{X^„U}.)VJ{Xi.,W,)dBl 

J[0,t]xNx[0,l]xC/2 I \ - ly r, rJI X ^y ,, rj _ j^ 

+ fiK^v)'^{ei^{x,.,Ur)<e<ei{Xr,Ur,v)} >Qidr, di, dO, dv), (3.2) 

where 6\, 9\ and 0\ have been defined previously. 

By Remark 1 1.11 the knowledge of {ft, f) for f £ Cq' is enough to characterized the 
finite measure- valued process v. 
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We introduce the canonical filtration 

Tt = cj{v^-, B^JGN*; Q([0,r]xA), A G V{N) ® B{[0,1] xU), r <t}, 
where ^([0, 1] x U) is the Borel a-field on [0, 1] x U. 
Lemma 3.1 The measure 

q{dt, di, de, dv) = CsJ2 hn<t<n+,}{^k + ^)Y1 ^Mi)dtde—^dv 

= Cs{Nt + l)f^6,{dz)dtde^^dv 

is the (predictable) compensator of the multivariate point process Q . 

Proof. For A; > 0, a regular version of the conditional law of {Tf^+i, Ik+i, ^fc+ii ^fe+i) with 
respect to cr{i^Q, {B^), j £ W, {Tp, Ip, Op, Vp), 1 <p <k} is given by the measure 

CsiN, + l)l|^,<,|e-^^^^(^^+i)(*-^^) Y:S,{d^)dtde^dv. 

The conclusion is thus a consequence of |21j Theorem 1.33 p. 136. □ 

Using Lemma 13.11 and Ito's formula, one can immediately show that any solution 
v of (|3.2() . such that i?(sup«j^(f(, 1)^) < +cxd, is a Markov process with infinitesimal 
generator L defined by ()2.11|1 . Moreover, we also deduce the following existence, moment 
and martingale properties. 

Proposition 3.2 1) Assum^e Hypotheses (H) and that E{{i'q,1)) < +oo. 

Then E{si\pf^j-{i't, 1)) < +oo for each T > and the process v defined by \3.°<3\) is 
well defined on M"*". 

2) If furthermore for some p > 1, E{{vq, 1)^) < +oo, then for each T > 

^(sup(i/i,l)P) < +00. 
t<T 

Proof. We firstly prove 2). For each integer fc, define Sk = inf{t > 0, (i^t,l) > A;}. A 
simple computation using 1)3. 2() . and dropping the non-positive death terms, gives 

/ f-thSk \ 

E{ sup {v,,IY)<e({uo,IY + C {l + {u,,lY)ds] 

sG[0,tAS'fe] \ Jo J 

11 



Gronwall's lemma implies that for any T > 0, there exists a constant C independent of k, 
such that £'(sup(groT/\5^i(fi, 1)''') < C. One easily deduces that 5"^ tends a.s. to infinity 
when k tends to infinity and next, Fatou's lemma yields E(sup)-iz\QT]{^t, ^)^) < +oo. 

Point 1) is a consequence of point 2). Indeed, one builds the solution {ut)t>o step by 
step. One only has to check that the sequence of jump instants {Tk)k goes to infinity a.s. 
as k tends to infinity. But this follows from E{sup^^rp{ut, 1)) < +oo . □ 

The following martingale properties are the key point to study large population ap- 
proximations. 

Theorem 3.3 Assume Hypotheses (H) and that for some p >2, E{{i^q, 1)^) < +oo. 

1) Then, for F and f £ Cq'° such that for all u £ M, |F/(i/)| + |LF/(i/)| < C(l + {ly, 1)p), 
the process 

Ffii^t) - Ffii^o) - f LFf{vs)ds 
Jo 

is a cddldg martingale starting from 0. It is in particular true for F{y) = yP~^ . 

2) The process Z^ defined for f £ C^ by 

^t = Wt,f) - (i^oJ) - / \m{x,u)Axf{x,u) + b{x,u).V^f{x,u) 

is a cddldg L'^ -martingale starting from with predictable quadratic variation 
{Zf)t= I [ {2m{x,u)\V,f\^ + {X{x,u)+fi{x,u,I^W*Us{x,u)))f^{x,u) 

Jo JxxU I 



+ f {x,v)M{x,u,v)dv >i's{dx,du)ds (3.4) 

Proof. Point 1) is immediate. For point 2), we first assume that E({i'q, 1)^) < +cxd. 
Applying point 1) with F{y) = y (or ()3.2() and Lemma l3.1|) leads to Z-' . Then one applies 
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1) again with F{y) = y^, and thus 

Wtjf-ii^ojf- / \2{m{x,u)A^f{x,u) + b{x,u).Vxf{x,u)){iys,f) 
Jo JXxU I 

+ 2m{x, n)| V,/|2 + \{x, n)(2(z.„ /)/(x, u) + f\x, u)) 

+ / (2/(x, t;)(i/„ /) + f{x, v))M{x, u, v)dv 
Ju 

+ ^{x,u,I W -k iys(x,u)){—2{iys, f)f{x,u) + f {x,u))>i's{dx,du)ds (3.5) 

is a cadlag martingale. In another hand, Ito's formula allows us to compute {ut, f)'^ from 
(|3.3|) : the process 

{^tJ?-{yoJ?-l I |2(m(x,u)A,/(x,n) + 6(x,n).V,/(x,n))^„/) 

+ 2(A(x, n) — ii{x, u, I^W * i^s{x, u))){us, f)f{x, u) 

2f{x,v){i's, f)M{x,u,v)dv >Vs{dx, du)ds — {Z-')t (3.6) 

lu J 

is a cadlag martingale. Comparing ()3.5|) and H3.6|) leads to p. 4(1 . The extension to the 
case where E{{vq, 1)^) < +cxd is straightforward, noticing that E{{Zf)t) < +oo. □ 



4 Large population approximation for a fixed interaction 
range 

We are now interested in deterministic approximations of the population point process 
when the size of the population increases. We assume in this section that the interaction 
range (5 > is fixed. 

Let us consider a sequence of initial measures {iyQ)NeN* belonging to Ai. For each 

A^ G N*, we keep all parameters {m,b,X,M) unchanged, except the competition kernel. 

We assume that for each A^, 

r 
I2n{x, u, r) = ^{x, u, —). (4.1) 

This assumption has a biological interpretation. In a case of fixed amount of available 
global resources, a large system of individuals may only exist if the biomass of each in- 
teracting individual scales as j^, which implies that the interaction effect between two 
individuals scales as -^ as well. The parameter N can also be interpreted as scaling the 
resources available, so that the renormalization of /i reflects the decrease of competition 
for resources. 
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We assume that the sequence ^ converges, as A^ tends to infinity. The size {i'q , 1) 
of the population is then of order N and will stay at this order (or at a smaller order) 
during finite time-intervals, since birth rates are bounded. Hence, our aim is to study the 
asymptotic behavior, as A^ tends to infinity, of the cadlag process 

Af = ^i;i(x;«, = >.", (4.2) 

taking values in M^ = {j^J^, v G M^. 

The process (A^)t>o is a Markov process with generator Ljv = ^n,\ + -^Af,2- An easy 
computation, for F £ C^(M) and / G Cq' , gives that 

i>7V,2i^/M = (^, m(.)A,,/ + 6(.).V,/)F'((z., /)) + ^, !^|V,/|2)F"((z., /)) (4.3) 
and (using (|4.1|) ) 

L^,iF;(z.) = iV / |A(x,n)(F((z.,/) + ^/(x,n)) - F{l,v,j))) 

+ ^x{x, u, I'W * u{x, u)){F{{u, f) - l/(x, n)) - F((z., /))) 

+ J^(F((z.,/) + l/(x,i;))-F((i.,/)))M(x,n,«)(it;|z.(dx,du) (4.4) 

We deduce from Theorem 13.31 the following martingale properties. 

Lemma 4.1 Let N > 1 be fixed and assume that for some p >2, E ((A^, l\ ) < cxd. For 
all Cq' -function f, the process 



zf'^ = (Af,/>-(A^,/>- / / \m{x,u)AJix,u)+bix,u).VJ{x,u) 

Jo JxxU I 

+ (A(x,u) - fj.{x,u, I^W * A^ {x,u))) f{x,u) + / f{x,v)M{x,u,v)dv >A^{dx,du)ds 

(4.5) 

is a cadlag Lp' martingale starting from with predictable quadratic variation 

(Z^'/)i = l / / {2m[x,u)\VJ\^ + {X{x,u) + ^l{x,u,I'W ^ A^ (x,u)))f\x,u) 
^^ Jo JxxU I 

+ / f'^{x,v)M{x,u,v)dv\A^{dx,du)ds (4.6) 

We assume 

Assumption (HI): 

1) The initial measures Aq converge in law and for the weak topology on Mf{X x U) to 
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some deterministic finite measure ^q E Mp{X x lA), and supjy E{{Aq , 1)^) < +cxd. 

2) All the parameters of the model are assumed to be continuous, either on X xU, or on 
X xUy.^. 

3) There exists a constant k^ such that 

Vx e A", u G Z^, ri,r2 G M, |/x(j;, u, ri) - /i(x,u, r2)| < A;^|ri - r2|. (4.7) 



By the law of large numbers, Assumption (Hl-1) is for example satisfied for Aq = 
Iv X^j=i "^(xSc/Mj with independent random variables (ATq, ^o){i<j<Ar} distributed follow- 
ing the law ^0 with finite 3rd-order moment. 

Let us recall that the parameters of diffusion, birth and mutation associated with A^ 
stay unchanged, whereas the parameter of selection /i^v is defined by (|4.1j) . 

Theorem 4.2 Assume Hypotheses (H) and (HI), and consider the sequence of processes 
A^ defined by \4-^ - Then for all T > 0, the sequence (A^) converges in law, in 
D([0, T], Mf{Xx lA)), to a deterministic continuous function ^ belonging to C([0, T], Mf{Xx 
U)). 

This measure-valued function ^ is the unique weak solution satisfying sup(grQTi(^j , 1) < 
+00 of the following nonlinear integro- differential equation. For all function f £ Cq' , 

iCf, f) = (eo, f)+ [ [ l^x, u)A,f{x, u) + 6(x, n).V,/(x, u) 
Jo JxxU I 



+ {\{x,u) — (i{x,u,I W -i^ ^g{x,u))^f{x,u) + I f{x,v)M{x,u,v)dv> S,s{dx,du)ds 

(4.8) 

Remark 4.3 Applying i4-^ to the constant function equal to 1, the positivity of fi and 
Hypotheses (H) gives {S,t, 1) < ('^O) 1) + C fQ{^s, l)ds. We conclude by GronwaWs lemma 
that any solution ^ of \4-^ i^ bounded on every finite time interval [0,T]; 

sup (ef,l) < (^o,l)e^^. 
te[o,T] 

As a first step in the proof of Theorem 4.2, we now give a mild formulation for solutions 
of 1)4. 8() . To this aim, and for each fixed trait u £ U, we denote by P" the semigroup 
of the diffusion process normally reflected at the boundary of X, with diffusion matrix 
m(-,u)Id and drift coefficient b{-,u). 

Lemma 4.4 Let us consider a solution ^ of i4-^ - Then, for each measurable and 
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bounded function (f defined on X xU, 

r-t 



Jo JXxU I 



+ / Pl'_,(f{x,v)M{x,u,v)dv^ Cs{dx,du)ds. (4.9) 

Proof. We may classically derive from ()4.8() a space-time weak equation for measurable 
functions iIjs{x, u) = ■(/'(s, x, u) which are of class C^'"^ on [0, t] x X, measurable and bounded 
in u and such that dnip = on [0, t] x dX x lA, given by 

{it^i^t) = {io^i^o) + / / \dsilJs{x,u) +m{x,u)l:^x'^s{x,u) +b{x,u).VxiJs{x,u) 



Jo JXxU L 
+ (A(x,ti) — /i(x,u, / W •k^g{x,u)))'ips{x,u) + I '4Js{x,v)M{x,u,v)dv\ ^g{dx,du)ds 

(4.10) 

Let us now consider a continuous function ip on X x U and fix a time t E [0,T]. Let us 
define for (s, x, u) G [0, t] x ^ x K, 

1/Js{x,u) =Pt_gip{x,u). 

Then f/; is solution of the boundary value problem 

dgipsixju) + m{x , u) Axipsix , u) + b{x,u).'Vx'4's{x,u) = on [0,T] x X x U 
dntpsix, u) = on [0, T]xdX xU 
iptix, u) = ip{x, u) on A^ X U. 

Equation (|4.1U|) applied to this function tp yields the evolution equation 

(ef,(^) = (eo,i't»+ / / \{Xix,u)-f,ix,u,l'W^Ciix,u)))P,^_Mx,u) 

Jo JXxV I 



+ / Pt_s(p{x,v)M{x,u,v)dv^ ^'^g{dx,du)ds. (4.11) 

Equation H4.11() is true for each continuous (and then bounded) function ip, and charac- 
terizes the finite measure ^ . Lemma 14.41 is proved. □ 

Proof, (of Theorem 13^. Let us fix T > 0. 

Let us firstly prove the uniqueness of solutions ^ of (|4.8jl . Using Remark 14.31 and Lemma 
14.41 we prove the uniqueness of bounded solutions of 1)4.91) . Let us consider two such 
solutions {£,t)t>o and {^t)t>o and compute the quantity |(^i — S,t,v)\, for each measurable 
and bounded function c/? such that ||(/?||oo < 1- 
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Using (|4.9|1 . we obtain for t <T 



(6-6,^)1 < 



XxU 



(A(x, u) - n{x, u, fW -k is{x, u))) Pt.s^{x, u) 

+ I Pl'_gip{x,v)M{x,u,v)dv> {^s{dx , du) — ^s{dx , du)) 



ds 



+ 



[fj,{x , u, I^W * ^s{x , u)) — n{x , u, I^W -k S,six , u))) P^_^g(f{x , u) S,s{dx,du)d. 



10 JXxU 

Now, using Hypotheses (H), applying Remark 14.31 to ^ and since ||(/?||oo < Ij there exists a 
positive constant Ci such that for all {x^u) ^ X xli and all < s < t < T, 

\\{x,u)P'^_,^{x,u) + P^_,^{x,v)M{x,u,v)dv\ < Ci, 
Ju 

\fiix,u,I^W*Ux,u))Pr^s^ix,u)\<fioil + \\I^W\\^{^s,l)) < Ci 



while thanks to (Hl-2), 

\n{x,u,I^W-k^s{x,u))- fi{x,u,I^W-kCsix,u))\ <k^\\I^W\\oo sup 1(6 - 6, V?)!, 



l¥'l|oo<l 



and then 



XxU 



(fi{x,u,I^W-k^s{x,u)) — fi{x,u,I^W-k£,s{x,u))) P^_g(p{x ,u)^s{dx, du) 



<C2 sup |(6-6,¥')l 

||¥'||oo<l 

where C2 is a positive constant. We deduce that there exists C > such that 

I {it -^t,^)\<C sup I {Cs -Cs,^) \ds 

Jo ||l(3||oo<l 

and by Gronwall's lemma, we conclude that for all t < T, supn 11 <]^ |(6 — 6,^)1 = 0. 
Thus, for all t <T, ^t = S,t and uniqueness holds. 

Let us next prove that for all T > 0, 



sup E sup (Af , 1>M < + 

A^eN* \[0,T] / 



00 



(4.12) 



Introducing S^ = inf{t > 0, (A^, l\ > A;} for k G N*, a simple computation using the 
specific form of L^^iFf and L]\f^2Ff with / = 1 and F{y) = y^ and dropping the negative 
death term yields 



vN i\3 



^N i\3 



E sup (Af , 1>^ < Ei{A^, ly) + CEl ((Af , 1) + (Af , l>^)d. 



,s<tAS^ 



tASf 



^N i\3n 



where C is a positive constant independent of k and N. Then Assumption (Hl-1) and 
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Gronwall's lemma imply that there exists a constant Ct independent of k and N such that 
E (swpg^rp^gN (A^, l) j < Ct- We deduce that the sequence {S^)k tends a.s. to infinity 
and finally obtain ()4.f2|) by Fatou's lemma. 

Using Remark 11.11 and following Roelly [37| , one observes that the sequence of laws 
Q^ of A^ is uniformly tight in V{^{[Q,T],Mf{X x U))), where Mp is endowed with 
the vague topology, as soon as for any function f £ Cq' , the sequence of the laws of 
the processes (A^,/) is tight in P(B([0,r],M)). Using Aldous' [I] and Rebolledo's [H] 
criteria, this tightness follows from 

sup ^(sup|(Af,/)|)<oo, (4.13) 

NeN* [0,T] 

and from the tightness of the laws of {{Z^'f)) and of the drift part of the semimartingales 

(A^,/). 

Clearly, since / is bounded, ()4.13() is a consequence of 1)4. 12() . Let us now consider stopping 
times {S, S') satisfying a.s. Q<S<S'<S + 5<T. Thanks to Doob's inequality. Lemma 
O and (gmi), we get 

E {{Z^'f)s> - {Z^^f)s) < E {c j'^' ((Af , 1) + (Af , 1)') ds^ < CS. 

Similar arguments prove that the expectation of the finite variation part of ('A^,,/') — 
(A;^,/) is bounded by C6. Finally it turns out that the sequence {Q'^)n is uniformly 
tight. 

Let us now denote by Q the limiting law in 'P{B){[0, T],Mp{X x U))) of a subsequence 
of Q^ , still denoted by Q'^ for simplicity. By construction, almost surely, 

sup sup KAf,/)-(Af_,/)|<l/iV. 
te[o,T]||/||oo<i 

We deduce immediately that each process A with law Q is a.s. strongly continuous. Let 
us finally prove that it is the unique solution of l|4.8|) . 

For t < T, / G Co'° and u G O{[0,T],Mf{X xU)), let us define 
^ij(^) = (^t'Z) - (^0,/) - / / \m{x,u)A^f{x,u) + b{x,u).V^f{x,u) 

^....- -::i,.„„../,.,.„„ 4,». 

We want to show that for any t <T, 

Ei\%j{A)\)=0, (4.14) 
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knowing from Lemma l4. II that 

Zf '^ = ^t,/(A^). (4.15) 

A fair computation using Lemma l4.H Hypotheses (H) and (HI), and (|4.12|) shows that 

E (|zf '^P) = E ((Z^'/),) < ^i^ (/* {l + (Af , 1>'} ds^ < %i (4.16) 

which goes to as A^ tends to infinity. On another hand, since A is a.s. strongly continuous, 
since / G Cq' and thanks to the assumption (H), the function ^tj is a.s. continuous at 
A. Furthermore, for any v G D([0,T], Mf(A' x U)), 

2 



\^tj{y)\ <a,/sup l + (l/s,l 

[0,t] ^ 

and ()4.12|1 implies that the sequence {^t,f{^^))N is uniformly integrable. Thus 

\.raE{\^tj{K'')\)=E{\^tj{K)\) (4.17) 

Combining (Hl-1), (gUSJ), (|^TT^ and (|ITT7|) . we conclude that (|m|) holds and that (gSJ 
is satisfied for any / G Cq' . 

Then A is uniquely identified to ^ , and the sequence (A^) converges to ^ in P([0, T], Mp{Xx U)), 
where Mp{X x U) is endowed with the vague topology. To extend this result to the weak 
topology, we use a criterion proved in j28| . Since the limiting process is continuous, it 
suffices to prove that the sequence ((A^,l))jv converges in law to (? ,1) in D([0, T],R). 
We may apply what has been done above with / = 1. Theorem 14.21 is proved. □ 

In the next section, we will be interested in the limit of small spatial interaction range 
6. An intermediate result consists in proving the existence of a density for each measure 
i^f , t > 0. We make the additional 
Assumption (H2): 

1) The diffusion coefficient m{x,u) is of class C^ in x and the second derivative of m (in 
x) is a-Holderian, uniformly in u, with a > 0. Moreover, m is assumed to he positive. 
Hence, since X xU is a compact set, there exists rn,* > such that for all (x, u) £ X x U, 

m{x,u) > m* > 0. 

2) The drift coefficient b{x, u) is of class C^ in x and the derivative of b (in x) is a- 
Holderian, uniformly in u, with a > 0. 

Assumptions (H) and (H2) and the smoothness of dX allow us to adapt Sato-Ueno 
j4Uj (Theorem 2.1 and Appendix) to obtain the following lemma. 
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Lemma 4.5 There exists a unique function pt{x, u, y) defined on M_|_ x X x U x JY , con- 
tinuous in (t, X, y) and which is a density function in y G X such that for each continuous 
function (p defined on X xU , each {x,u) G X xU, 



P^ip{x,u)= / pt{x,u,y)ip{y,u)dy (4.18) 

Jx 

Let us now prove the propagation in time of the absolute continuity property of the 
measure- valued solution ^ . 

Theorem 4.6 Assume (H), (HI) and (H2) and that S,o{dx,du) = gQ{x,u)dxdu. Then for 
each time t, the measure ^f has a density g G L°°([0, T], L^) with respect to the Lebesgue 
measure on X xU. Moreover, for each t and u, the function gi{.,u) is continuous on X . 

Proof. Let us come back to the equation (|4.9|1 satisfied by ^ . 

Using basic results on linear parabolic equations, we construct by induction a sequence 
of functions ((7n)n satisfying in a weak sense 

a,5r^(x,n) = A,.(m(x,n)5r+'(x,n))-V,(6(x,n)<7r+'(x,u)) 

+ \{x, u)g'^{x, u)+ g'^ix, v)M{x, u, v)dv - fi{x, u, I^W -k g'^{x, u))g'^+^ (x, u) 
Ju 

5o+^(x,n) = gQ{x,u) 

V^5"+^(t,x,u).n(x) = yit,x,u)eR+xdX xK. (4.19) 

Thanks to the nonnegativity of go, fi, \ and M, and applying the maximum principle, we 
can show that the functions gn are nonnegative (see [H]). By symmetry of M, Equation 
()4.19() (understood in the weak sense) means that for all Cq ' -function (p from X xU into 

M, 

{9t^^,f) = (go,^) + / \{m{x,u)A^(p{x,u) + b{x,u).V^(p{x,u)]g'^~^^{x,u) 
Jo JXxU I V / 



+ I A(x,u)c/?(x,n) + / 99(x,u)M(x,u, f)(it; ) 5'"(x,n) 

- ^jl{x,u,I^W i^ g'^{x,u))ip{x,u)g'^^^{x,u) \dxduds. (4.20) 
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The associated mild equation writes as before: for each continuous function 99, 
{9t~^^i'^)= I ( _Ptix,u,y)(p{y,u)dy]go{x,u)dxdu 



Jo JXxU 



A(x,u)( / pt_s{x,u,y)ip{y,u)dy 

X 



+ / / pt-s{x,v,y)Lp{y,v)dy]M{x,u,v)dv 
Ju \Jx 



9six,u) 



-fi{x,u,I W-kg'^{x,u))l _pt-s{x,u,y)ip{y,u)dy\g^~^^{x,u)>dxduds. 

(4.21) 

Hypotheses on the coefficients allow us to apply Fubini's theorem and to obtain that 
for each (y, n) £ X x U, 

Jx 

+ / / \iKx,u)Pt-six,u,y)g^{x,u) + / pt-six,u,y)g]}{x,v)M{x,u,v)dv 
Jo Jx I Ju 

- ^{x, u, I^W * g^{x, u))pt-s{x, u, y)g^^^{x, u) \dxds. (4.22) 

Then, thanks to the nonnegativity of 5""*"^, we get 

+ / p(_s(j;, u, y)(7"(x,t;)M(x,M, f)^^ >(ixds, (4.23) 

and deduce easily, integrating over y S ^, using Fubini's Theorem, the symmetry of M 
and Gronwall's Lemma that there exists a constant C independent of b such that 

supsup||<7ri|i<ll5o||ie^^. (4.24) 

n&it<T 

Let us now show the convergence of the sequence g^ in L°°([0, T], L^) to a function g . 
A straightforward computation using (|4.22(1 . Hypotheses (H), (HI) and (H2), and similar 
arguments as above yields 

ft 

n+l nil ^ r^ I I on-nll^'^ + 1 ^'^11 _LoiTi^ll^" ^"~1| 



sup hT' - <?"lli <C sup llsr' - glWi + sup \Wu - 9u~'\\i]ds 

s<t Jo \ u<.s u<s J 

where C is a positive constant. Thanks to Gronwall's Lemma, we deduce that for each 
T > 0, each t < T and each n, sup^<t ||g^+^ - y"!]! < C'/q sup„<^ \\g'^ - g'^'^Wids. 
Picard's Lemma allows us to conclude that for any T > 0, X^^g^ suptgro,r] ll^"^ ~5rili < 
+CXD, and the sequence {g"')n converges in L°°([0, T], L^) to a function g . We deduce 
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from !(T?I^ that 



sup sup Ibf 111 < Ibollie^^- (4.25) 

5>0 t<T 



Moreover, the function g is solution of H4.8() . and thus, the uniqueness result proved 
in Theoreni l4. 21 implies that £, {dx, du) = g {x, u)dxdu. Then, the measure ^ is absolutely 
continuous with respect to the Lebesgue measure, and the density g is weak solution of 
the nonlocal nonlinear partial differential equation 

dtg^ = A,^{m{x,u)gf{x,u)) - V,j;{b{x,u)gf{x,u)) 

+ {X{x,u) - fi{x,uJ^W -k gf{x,u)))gf{x,u) + / gf{x,v)M{x,u,v)dv ; 

Ju 

9o(.x,u) = go{x,u) ; 

V^/(t, X, u).n{x) = V(t, X, u) £R+xdX X U. (4.26) 

Lemma 14.41 implies that g is also solution of the mild equation 

Jx 

+ UHx,u)-n{x,u,I^W-kg^^{x,u)))pt-six,u,y)gl{x,u) 

+ / pt-s{x,u,y)gg{x,v)M{x,u,v)dv>dxds. (4.27) 

Using (|4.27|) . the continuity of y i-^ g^{y,u) follows immediately from the continuity of 
{t,x,y) 1-^ pt{x,u,y), the nonnegativity and boundedness of g and the boundedness of 
birth parameters. □ 



5 Convergence of the number density when the interaction 
range decreases 

Our aim in this section is to prove that under suitable assumptions, the sequence {g ) 
converges, when 5 tends to 0, to a function g G L°°([0, T], L^) with initial condition gQ, 
which is weak solution of the locally nonlinear partial differential equation 

dtgt{x,u) = A^{m{x,u)gt{x,u)) -Vx{b{x,u)gt{x,u)) 

+ [X{x,u)-fi{x,u,pg{t,x,u)))gt{x,u)+ / gt{x,v)M{x,u,v)dv ; 

Ju 

Vxg{t, X, u).n{x) = V(t, x,u) £R+x dX xU (5.1) 
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where pg describes the (local) interaction in x, defined for {x^u) ^ X xU by 



Pg{t,x,u)= / W{u-v)gt{x,v)dv. 
Ju 

In order to control the terms I W * g uniformly in 5 in the nonlinear term of ()4.26|) . we 
need L°°-estimates on g and we make the following initial data assumption: 

(H3) The initial density gQ{x,u) is bounded on X xU. 

Proposition 5.1 Assume (H), (HI), (H2), (H3). Then there exists a positive constant 
Ct, such that 

sup sup lls-f lloo < C't||5o||oo- (5.2) 

s>o te[o,T] 

Proof. Let us again consider the sequence {g^)n approximating g introduced in the proof 
of Theorem 14.61 The maximum principle implies that 

SUpSUp||5"||oo < (^11501100, 
n t<T 

where C > is a constant only depending on T, A* and M* (and independent of 5). 
This property propagates taking the limit in n, and 1)5. 2|) is proved. (For details on the 
maximum principle, see [H].) □ 

Let us now prove the following convergence theorem: 

Theorem 5.2 Assume hypotheses (H), (HI), (H2), (H3). Assume that the measure 
I {y)dy weakly converges to the Dirac measure 6o as 6 tends to 0. (To fix ideas we may 
assume that I (x) = CsliixiKS}-) Then the sequence {g )s>o converges in L°°{[0,T],L^) 
as 6 tends to 0, to the unique function g G L°°([0, T],L^ Ci L°°{X x lA)) satisfying for each 
y,u & X X U the evolution equation 

9tiy,u)= / ptix,u,y)go{x,u)dx 
Jx 

+ UHx,u)-lJ'{x,u,pg{s,x,u)))pt-six,u,y)gsix,u) 



+ / pt-s{x,u,y)gs{x,v)M{x,u,v)dv>dxds. (5.3) 

Moreover, for each t and u, the function gt{.,u) is continuous on X. 

Proof. One can easily prove the existence and uniqueness of the integrable and bounded 
function g solution of 1)5. 3|) by adapting the proofs of Theorem 14.61 and Proposition 15.11 
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replacing ij{x,u,I^W -k g) by fi{x,u,pg). The continuity oi y —>■ gt{y,u) is obtained as in 
the proof of Theorem 14.61 and we can show as in the proof of Proposition 15.11 that 



sup \\gt\\oo < CtWoqI 
te[o,T] 



(5.4) 



Let us write 



9tiy,u) - gtiy,u) 



Jx 



(A(x, u)pt-s{x, u, y) {gi{x, u) - gs{x, u)) 



+ / pt-s{x,u,y){gs{x,v) - gs{x,v))M{x,u,v)dv 
Ju 

fi{x,u,I^W -k gl{x,u))gl{x,u) - iJ.{x,u,pg{s,x,u))gs{x,u) pt-s{x,u,y) \dxds 

(5.5) 

Using (|5.2|) and (|5.4|) . the unique term which deserves attention is the term /i(x, u, Pg{x^ u))- 
IJ,{x,u,I^W -kg^{x,u)). By (|I?7)) . we have 

\p{x,u, pg{t, x,u)) — p{x,u, I W -k g^ (x, u))\dx 



X 

< k,. 



X 



X 



W{u — v)gt{x,v)dv — / / {x — z)W{u — v)gf{z,v)dzdv 
u JxxU 



dx 



W{u — v)gt{x,v)dv — / / {x — z)W{u — v)gt{z,v)dzdv 
u JxxU 



+ / I%x-z)W{u-v) gt{z,v)-gtiz,v) 
iXxU 



dzdv j dx 



Let us fix our attention on the first term in the last right inequality, that we will call 
A^{t, u). Since I^{y)dy weakly converges to 6o, J-^ I^{x — z)gt{z, v)dz converges to gt{x, v) 
as 5 goes to 0. Because of (|5.4() . this convergence holds in a bounded pointwise sense 
with respect to t < T, x G A" and v £ U. Then Lebesgue's theorem implies that A^^x '■= 
jjj Jq As(t, u)dtdu tends to as 5 tends to 0. 

Now, integrating 1)5. 5() with respect to dy du, a straightforward computation yields 
sup \\gi - gs\\i< CtAs^t + 0^ sup \\gi - gu\\ids. 

s<t Jo u<s 



We conclude using Gronwall's lemma. 
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The zero interaction range equation 1)5. 1|) has been numerically studied in Prevost |34j . 
A lot of simulations based on finite element schemes are given, studying the simultaneous 
effects of the diffusion, mutation and selection on the invasion of the domain by the 
population. The simulations show that the coefficient which seems to affect the most 
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the invasion aptitude is the mutation size coefficient. However, they restrict to local 
interactions. 

In the next section, we wish additionally to illustrate, by simulations of the stochastic 
discrete model, the effect of the spatial interaction range on the interplay between invasion 
and evolution, and the emergence of spatial and phenotypic diversity (clustering and 
polymorphism). Our simulations focus on the qualitative differences between local and 
nonlocal interactions. 

6 Simulations 

We will give in this section simulations of several biologically realistic examples, based 
on the algorithm of Section |H| The Euler scheme to simulate reflected diffusions will be 
detailed in Section [6. 11 as well as some simpliflcations in the algorithm of Section |31 in the 
case of linear death rates. 

Next, we will give simulations of three biologically relevant examples. First (Sec- 
tion |ni2)i we show that, when migrations and mutations are not too strong, a large in- 
teraction range induces a spatial organization of the population as a flnite set of isolated 
clusters. Conversely, for sufficiently small interaction range, the clustering phenomenon 
is no more observed. Second (Section 16. 3() . we propose another example where a similar 
phase transition occurs for spatial clustering and in which the critical interaction range 
can be identifled. In our last example (Section 16.4(1 . we investigate a model describing 
the invasion of a species with evolving dispersal speed (the trait is proportional to the 
migration speed, as in [HI). 

6.1 Euler scheme and algorithm for logistic interaction 

As mentionned in Section |HJ the reflected diffusion of our particles can be simulated with 
an Euler scheme. We will assume in this subsection and in the following examples that 
A" = (a,/5) C M and we will use the scheme of Lepingle j21] (see also HZ]). Fix x G [a, 13] 
and u £U. On any time interval where its trait is constant, an indivual at {x,u) moves 
according to the reflected diffusion 

Xt = x+ / ^/2m{Xs,u)dBs + I b{Xs,u)ds - h , 
Jo Jo 

\k\t = I '^{x,&{a,i3}}d\k\s ; kt= (l{X,=/3} - l{x,=a}) d\k\s, (6.1) 

JO Jo 

where i? is a one-dimensional Brownian motion. 
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If m and b are Lipschitz with respect to the first variable, then one can simulate this 
diffusion on [0, T] as follows. Fix h > 0, a and /? such that a < a < JS < (3, and let n be 
the first integer greater than T/h. For p € {0, 1, . . . , n — 1} and ph < t < {p + l)h, let 

Xo = X, 



Xt = max[a, min[/3, Xph + b{Xph, u){t - ph) + y2m{Xph, u){Bt - Bph) 

+ l{Xp„<s} max(0, A^ - {Xph - a)) 
- l|x,,>/3} max(0, i?f + {Xph - /?))]], 



where 



A^. = sup <^ -b{Xph, u){s - ph) - J2m{Xph,u){Bs - Bph) \ , 

ph<s<t I J 

5f= sup {b{Xph,u){s - ph) + J2m{Xph,u){Bs - Bph)] . 

ph<s<t L J 

Then, there exists a constant C independent of h such that for any function / on [a, /?] 
with finite variation bounded by 1, supQ<j<j^ \E{f{Xt) — f{Xt))\ < C\fh. 

In each step of this scheme, one has to simulate simultaneously Bt and St := sup^<j(ai?s+ 
bs) for fixed constants a, b. This can be done as follows (Shepp [41j). Let C/ be a Gaussian 
centered random variable with variance t, and let V be an exponential random variable 
with parameter l/2t independent of U . Put 
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y = - UC/ + 6t + ^d^V + (aC/ + 6t)2 

Then the vectors {Bt,St) and (t/, i^) have the same distribution. 

Note that this scheme can be easily generalized to state spaces of the form (ai,/3i) x 
... X {ad, Pd) C M*^, as explained in ^^. 

Next, we want to study a particular case in which we can considerably reduce the 
complexity of the algorithm. In Section |31 one needs to compute I W * v{x,u) at some 
point {x,u) £ X X U at each time step, which involves a sum over all individuals in the 
population. In the case of logistic competition (linar death rate) where 

fi{x,u,r) = fio{x,u) + tii{x,u)r, (6.2) 

one can use the following algorithm. 

Fix a constant Cs in a similar way as in ()2.6|) such that pq{x,u) + X* + ||A/*||i < Cs and 
/ii(x,u)||/ I^lloo < Cs- Take the Brownian motions {B^' )j^keN and the random variables 
(^fc)fcgN) {Vk)keN and {Tf:)keN as in Section |21 Set Tq = and Nq = N (the initial number 
of individuals). Assume that {Tk-i, Nk~i, Xt^_-^, Uti^_-^) are given. N^-i is the number of 
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individuals at time Tfc„i. At this time, their positions and traits are the coordinates of 
the vectors Xt^_^ = {X^_ )i<i<Nt,_i and Ut^-^ = {U^,_ )i<i<N^^i- The two first steps of 
the algorithm are the same: the new time step is given by T^ = Tf^-i + T/fc/A'^fc_i(A^^._i + 1) 
and the motion of each particle is governed by the SDE with normal reflection H3.1() . 

The third step deals with the different events that may happen at time T^. Choose 
at random one individual Ik = i uniformly among the N^-i individuals living during the 
time interval [Tfc_i,Tfc). Its position and trait are {X^ ,U^ ). The event occurring at 
time Tfc is decided by comparing 9k with constants related to the rate of each kind of 
event. The only difference with the algorithm of Section |31 is in the first sub-step, that has 
to be divided in two steps as follows: 

• If < ^fc < j^^ =: 91{Xt„Ut,_,), then let j G {1, . . . , iV^.i} be such that 

J-1 < R, ^ J Tf a, _ J-1 ^ ^^ ^fc ^k-l' ^ ^fc T^' ^ ^fc-1 T^-l' 

then the individual i dies from competition with individual j and N^ = N^-i — 1- 
Otherwise, nothing happens and N^ = Nk-i- 

. If ^^(Xt,, t/T,_J <0k< ei{XT,,UT,_,) + cs$jl'i] =■■ OHXt^Ut,^,), then the 
individual i dies by natural death and Nk = Nk^i — 1. 

The three other sub-steps are the same. 

The main difference with the algorithm of Section El is that we no longer have to 
compute ^j=i^ ^^(-^Tk ~ -^rJ^i^T^^i ~ ^T^^i) ™ ^^^ ^^^^ sub-step, but it suffices to 
compute I {X^ — X!^ )W{U^ _ ~^t ) ^^^ chosen i and j. Moreover, we do not need 
to compute the position of each individual in the population at each time step. The third 
step above only needs to compute the position of at most two particles at time T^ (the 
particles numbered i and j). 

We assume in the following examples a logistic competition of the form ()6.2|) and a 
physical space of the form X = (a,/3). Our simulations are realized with the previous 
algorithm. 

6.2 Example 1. Spatial clustering 

We consider here a set of parameters similar to the one of jlUj and j34j . in which, for 
each spatial position x, the growth rate is maximal for the trait value u = x. This can 
represent the effect of a gradual spatial distribution of different resources, involving a 
gradual distribution of traits. For example, for some bird species, a linearly spatially 
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varying seed size determines a linear variation of the beak sizes (Grant and Grant ;18j). 

A" = (0,1), U=[0,1], m{x,u) = m, b{x,u) = 0, 

X{x, u) =2 — 20(x — u) if Ix — n| < 1/vlO; otherwise, 

//Ar(x,u,r) = 1 + — , /(y) = C5l{|y|<5}, W{v) = l. 

Moreover, M{x,u,v) = 0.1 x ks{u,v) where 0.1 is the mutation rate and ks(u,v) is the 
probabihty density of a Gaussian random variable with mean u and standard deviation 
■s conditioned on staying in U = [0,1]. Therefore, we have four free parameters in this 
model, m, 6, s and the population size N. The initial population in our simulations is 
composed of N individuals at (0.5, 0.5). 

The simulations of this model show, as in [23 ? that the invasion of space occurs along 
the diagonal x = u, and, as in ^Uj, that speciation (stable coexistence of several sub- 
populations with different typical traits) may occur in this model, accompanied with 
a spatial specialization. Several different population clusters may coexist at different 
position, with trait values located around the corresponding optimal traits. We have 
investigated in our simulations the effect of the different parameters on the clustering and 
polymorphism of the population. We give pictures of the seemingly stable state of the 
population (Fig.^. Our first general observation is that the clusters are more concentrated 
at the boundary of the domain. Indeed, the reflected diffusion governing the motion of 
individuals is not isotropic close to this boundary, so that the population density is bigger. 

On the one hand, we have investigated the effect of the interaction range 6. The 
main result is that the existence, the number of clusters and the distance between clusters 
strongly depend on the interaction range 6. As shown in Fig.^a-c), the number of clusters 
decreases with 6 and the distance between two population peaks is roughly 2S, which is 
exactly the width of the interaction interval. The emergence of population clusters is 
mainly a consequence of local births (HS^). Indeed, since the progeny of an individual 
is close to its original location, each individual's progeny can create a colony with stable 
position on short timescales. Once an individual is at a distance greater than 6 from the 
main part of the population, it experiences very little competition and it can create a new 
colony. When several colonies appear, they organize in a way to minimize the competition 
between them and to maximize the growth rate. 

If 6 is sufficiently small, we observe a flat distribution of the population (Fig.^c)), and 
thus a qualitative difference with respect to cases (a-b). As an explanation, decreasing 
6 increases the number of clusters and the width of a cluster increases with the speed of 
dispersal m and the range of mutation s. Then no distinct colony can be observed for 
sufficiently small 6, and fixed m and s. 
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(a) iV = 3000, s = 0.01, m = 0.01, S = 0.3. 



(b) iV = 3000, s = 0.01, m = 0.01, 5 = 0.1. 
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(c) iV = 3000, s = 0.01, m = 0.01, 5 = 0.03. 



(d) iV = 50, s = 0.01, m = 0.01, 5 = 0.1. 



"SSmatr160.txt" matrix 



■'SSmatr100.txt" matrix 





(e) N = 3000, s = 0.003, m = 0.03, 5 = 0.1. 



(f) iV = 3000, s = 0.03, m = 0.03, 5 = 0.1. 



Figure 1: Simulations of example 1 for various parameters. All of them are taken at time 
4000. 
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We also investigated the effect of the population size N. It appears that this parameter 
has very little qualitative effect on the clustering of the population, except for small 
N (Fig. n^d)), where the width of each clusters is reduced, and we observe much more 
fluctuations in the population distribution. However, we still can observe a similar pattern 
of population clusters than in Fig. ^b). 

On the other hand, we also studied the effect of the diffusion coefficient m and the 
mutation range s. Comparing Fig. ^b) and (f), we observe that too large s and m induce 
the same flat distribution as for small 6. This confirms that the clusters pattern depends 
mainly on the balance between m and s, and 5. In Fig. ^f), quick movements mix the 
population so that no spatial structure can appear. 

Finally, we also studied the relative effect of s and m in the appearance of spatial or 
phenotypic structure. As shown in Fig. ^e), small s can induce a differentiation over the 
phenotype space U even when m is large enough to have a fiat distribution over space X 
(compare with Fig.^f)). Fig. ^e) can be seen as an intermediate state between Fig. ^a) 
and Fig. ^^f). When m is reduced instead of s, a reversed pattern can be observed. 

6.3 Example 2. The role of spatial competition for clustering 

As we have seen above, the balance between the spatial competition range 6 and the 
diffusion parameters s and m has an important effect on the clustering of the population. 
Here we want to address the balance between the range of competition and the growth 
rate. For this purpose, we consider the following model, inspired by the adaptive dynamics 
model of [7]: 

X = {-1,1), U=[0,2], m{x,u) = m, b{x,u) = 0, 

\{x, u) = exp(-x^/2/9^), fiNix, u, r) = 1 + — , 

l\y) = Cs exp(-yV25'), W{v) = exp(-t;V0.02). 

and the same mutation kernel as above. This example has five free parameters m,6,s, N, 
and p, which represents the width of the space region with significant growth rate (namely, 
a parameter describing the width of the space region with high concentration of resources) . 
The initial population in our simulations is composed of N individuals at (0, 1). Observe 
that in this example, the trait has no effect on the growth rate, so that the trait structure 
is neutral (all individual's parameters are equal, independent of the trait, in absence of 
interaction) . 

Remark that if we consider the space A:" as a trait space, this model is similar to 
the one of 0. In particular, the biological theory of adaptive dynamics (^HJ) suggests 
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that evolutionary branching, i.e. the spht of the population into two sub-populations with 
different traits stably coexisting, translating in our model into spatial clustering, occurs 
when the range of interaction (6 in our case) is smaller than the range of the growth rate 
(p in our case). This is illustrated by Fig. Efa) and (b), where, m (a), 5 < p and the 
population stabilizes around two distinct positions (branching occurs) and in (b), 6 > p 
and the population stabilizes around position (there is no branching). 

Figures (c) and (d) prove that other phase transitions occur for smaller S, leading to 
the coexistence of three clusters or more. As in example 1, we notice in Figure (d) that 
very small 6 leads to a distribution without distinct clusters. 

It is possible to add some phenotypic structure to this example by assuming that the 
growth rate depends on the trait u, in a way such that spatial branching occurs for some 
traits but not for others, according to the above branching criterion. We take the same 
parameters, except for the birth rate, which has the following form. 

X{x, u) = exp(-xV2(n + 0.1)). 

The parameter p is then replaced by \/u + 0.1, so that branching occurs if \/u + 0.1 > 6. 

This is what happens actually for small times (Fig. Ol^a)), but after a longer time 
(Fig. |21^b,c)), the two clusters spread over the trait space because of the mutations. Even- 
tually, if we let time go on, we actually observe the appearance and the spread of three 
spatial clusters, colonizing all the trait space (Fig. I^Kd)). 

6.4 Example 3. Invasion and evolution of migration speed 

Here, we investigate a model describing the invasion of a species with evolving dispersal 
speed (as in [H]). This can model phenomena such as the invasion of Australia by cane 
toads, for which an adaptation to high invasion speeds has been recently detected (Phillips 
et al. |33j). The parameters are as follows. 

X = {-1,1), Z^=[0,3], m{x,u)=m{u + 0.1), b{x,u)=0, 

r 
X{x, u) = 1, pn{x, u,r) = 1 + —, 

I^y) = Ql||^|<5|, W{v) = exp(-107;2). 

and the same mutation kernel as above. Here we study invasion into an homogeneous 
space (A is constant). The diffusion rate m is proportional (up to a constant) to the 
trait u. Thus, individuals with large u move fast. The trait u can be a morphological 
trait responsible for the speed of dispersal (e.g. the length of legs for toads, [SSI)- Space 
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(a) 5 = 0.9, p= 1. 



(h) S = 1.1, p ^ 1. 
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"SSmatr039.txt" matrix 





(c) <5 = 0.5, p 



{d)S = 0.1,p- 



Figure 2: Simulations of example 2. Neutral case. All the figures are taken at time 5000, 
except the last one, taken at time fOOOO (more time is needed to fill the whole space). In 
ah the simulations, N = 1000, s = 0.003 and m = 0.003. 
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(a) t = 10000. 



(b) t = 20000. 
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(c) t = 45000. 



(d) t = 80000. 



Figure 3: Simulations of example 2. Trait-dependent case. In this simulation, A^ 
1000, s = 0.003, m = 0.003 and 6 = 1. 
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competition occurs between individuals within a distance 5, and the kernel W models 
competition between close traits. This example has four free parameters, the diffusion 
coefficient m, the interaction range 6, the standard deviation of mutations s and the 
population size A^. 

In Fig. |3] and El we present two extreme cases with respect to the initial trait distri- 
bution, but with identical parameters. In the first one, all individuals are at (physical) 
position 0, and with traits regularly distributed in U = [0,3]. In the second one, all 
individuals are initially located at a single point (0,0). 

In both figures, we observe a triangular invasion pattern indicating that the invasion 
front is composed of faster individuals. In Fig. |S1 we also observe a simultaneous invasion 
in physical and trait spaces, and a slower spread of the population. This can be explained 
by the fact that the population, initially composed of slow individuals, has first to colonize 
the trait space before invading the physical space. Because of the progressive appearance 
of larger traits, the invasion speed increases over time (compare the different time values 
in Fig. EJ. 

When parameters vary, the simulations of this microscopic model can show different 
ways of colonization. As an illustration, we give an example (Fig.|^ where the interaction 
range 6 is bigger. The parameters A^ and m are chosen such that two clusters appear for 
large traits and spread over the trait space in a short time. The initial condition is the 
same as in Fig. 0] We can observe two branches linking the initial cluster with the two 
extreme valued clusters (Fig. |Hl^c,d)). 
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Figure 4: Simulations of example 3. The parameters are A'^ = 100, s = 0.03, m = 0.003 
and 5 = 0.1. The initial condition is composed of N individuals located at and with 
trait values 3i/N for 1 < i < N. 
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Figure 5: Simulations of example 3. The parameters are A^ = 100, s = 0.03, m = 0.003 
and 5 = 0.1. The initial condition is composed of N individuals located at (0,0). 
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